#install.packages("dslabs","wooldridge","tidyverse") #Si ya están instalados, comenten con '#' para desactivar esta líneaCACCS: Taller de RStudio - 1ª parte
1 Introducción
Introducción a la computación estadística en R
Este taller está basado en notas de Paul Thibodeau y revisiones de profesores del departamento de psicología de BYU (vea esa versión aquí). Adaptado al español y usando ejemplos del libro de Rafael Irizarry disponible aquí. Ha sido revisado de la primeras versión de este taller, brindadas en el Centro Académico de Cómputos de Ciencias Sociales en los pasados semestres.
Si aún no has instalado R, está disponible aquí. Acto seguido, baja RStudio. Puedes también ir a la nube en Posit Cloud.
2 Sobre aprender R
R es un lenguaje creado por estadísticos como ambiente interactivo para análisis de datos. En R pueden guardar su trabajo como una secuencia de comandos, conocida como un script, que se pueden ejecutar fácilmente en cualquier momento, con portabilidad. Estos scripts sirven como un registro del análisis que realizaron, una característica clave que facilita el trabajo reproducible. Si bien es un programa poderoso y flexible, es ciertamente más complicado que programados como el que puedan encontrar en SPSS o STATA, donde pueden señalar con el ratón una opción y ejecutarla. Por otro lado, R es gratuito y de código abierto. En adición es modular y las funcionalidades añadidas complementariamente por terceros también son gratuitos, incluyendo acceso temprano a los métodos y herramientas más recientes que se desarrollan para una amplia variedad de disciplinas, incluyendo la economía, ciencias políticas, sociología, planificación, ecología, geografía o la biología molecular, entre otros.
Mi idea tras esta secuencia de talleres sobre R es brindar una introducción para que puedan interactuar con R tal que los sentimientos de frustración o ansiedad que puedan surgir al trabajar con este lenguaje al usarlo en variedad de circunstancias. Este taller presume en general ciertos conocimientos básicos en estadística, pero abundará en cada parte que se trabaje aunque sea algo como la motivación para las funciones que ejecutaremos. Si no tienen experiencia en programación computacional, es posible que partes de este taller resulten confusas o no tengan mucho sentido. Recomiendo que aún así traten de entender lo que puedan, y luego al volver a este taller periódicamente al repasar el contenido, vayan comprendiendo (o trayendo preguntas que no se resuelvan) el material.
Es importante que sepan que la habilidad más útil de antemano es la de saber adónde ir cuando se atasquen. R ayuda bastante en esto. La consola de R, al escribir y ejecutar comandos, brinda retroalimentación inmediata. Hagan uso de esto de manera liberal mientras trabajan el taller, haciendo pequeñas modificaciones a los ejemplos que brinde, hasta que sientan que entienden qué está ocurriendo.
De una vez, cabe señalar que la consola brinda acceso a la función de buscar ayuda interna. Casi todas las funciones (y muchos de los datos pre-almacenados) tienen archivos que les acompañan describiendo qué son las funciones, y cómo usarlas. Pueden acceder a esta información escribiendo help(función) (o ?función), sustituyendo «función» por el nombre de la función que quieras conocer. Es importante leer estos archivos detenidamente la primera vez que te encuentres con una función, pero también (posiblemente más) importante es consultarlos con frecuencia. Si tienes una idea de qué quieres hacer, pero no sabes ni recuerdas la función exacta para esto, puedes buscar en los archivos por el término usando doble signo de interrogación (??regression).
Hay mucha información adicional en internet. Es difícil buscar en Google por R, ya que es una letra, pero sí se puede encontrar mucho sobre funcionalidades. Para búsqueda más certera se puede ir a StackOverflow, y buscar el tag de R ([R]) junto a tu pregunta o error. Relacionado está el StackExchange de estadísticas, que es como Stack Overflow pero enfocado en la madeja estadística.
Finalmente, antes de empezar, quiero hablar sobre errores. Tanto novatos como expertos encontrarán errores en su código de R. De suceder, al ir ejecutando un script o trabajo, cesará el proceso y saldrá impreso un mensaje de error en la consola. Esto puede ser frustrante, en especial al estar iniciando el aprendizaje, pues el error ocurre a menudo muy adentro de las especificidades de una función, y el mensaje de error no guarda relación con lo que el usuario quería hacer. Un error a evitar al empezar a aprender R es que el error es un disparate incomprensible, y resignarse a la frustración y desánimo. Resistan ese impulso; si bien el mensaje no será de inicio informativo, está diseñado para transmitir cierta información con claridad, y entender esa información es clave para resolver el problema (o cambiar de estrategias).
Esta es una secuencia de dos talleres, y en esta lección arrancaremos con lo básico de R, RStudio y el tidyverse, y tenemos la meta hoy de que al culminar las primeras dos horas de esta secuencia podamos:
Entender la lógica de RStudio.
Importar, crear, leer y manipular datos y objetos.
Inicios de la visualización y análisis descriptivo.
En el siguiente taller expandiremos sobre esta base, y profundizaremos en visualización avanzada, así como en la inferencia estadística, con introducciones en R sobre modelos lineales, jerárquicos y longitudinales, así como de los distintos métodos de verificación de presunciones de modelo. Si bien hay mucho más que se puede hacer con R (e.g. machine learning, extracción de datos de la web, procesamiento de cadenas y minería de textos, análisis de redes, presentaciones estilo powerpoint, interacción con Git, trabajos en Unix, etc.), la idea es iniciar la travesía en R, y señalar a dónde pueden buscar ahondar.
2.1 Usando RStudio (y presentando por encima, Quarto)
¿Cómo difieren R y RStudio? RStudio es una interfaz que yace encima del programa base R. Es además mucho más amigable y llevadera que R.
R, por su cuenta es una consola sencilla donde teclean y corren código.
RStudio tiene la Consola como panel, pero tiene otros paneles en adición que nos ayudan a ubicarnos, como Environment (Ambiente) - para que vean los objetos (conjuntos de datos y variables, o funciones) que hayan guardado o creado a través de la sesión - el de Output (Salida) - donde pueden ver sus archivos en el directorio, sus gráficos, acceder los paquetes y buscar ayuda - y Source (Fuente), donde pueden crear y editar tanto archivos o scripts de R, como de Markdown, Quarto y Shiny (para dashboards), entre otros formatos.
Estos paneles (y otras configuraciones) son editables, como el caso del R que manejo en mi Mac:
A grosso modo, RMarkdown/Quarto es un tipo de documento que permite la escritura de texto regular junto a secciones o pedazos de R. Esto es útil para explicar por ejemplo tu código, seguido por el código mismo para uso regular. También es útil porque Markdown/Quarto se puede exportar como una página html así como un PDF, o Word lo cual es útil para compartirlo fuera del esquema de R. Si les interesara más sobre RMarkdown, visiten esta página o esta para Quarto aquí.
3 Interacción con la consola de R
En la forma básica, R es una calculadora para cómputos básicos. Escribe en R (o en un trozo o chunk en RMarkdown o Quarto) y el resultado saldrá en la consola al ejecutar o correrlo. Un chunk se designa en Markdown con triple ` al inicio y fin, así como {r título de sección}. Ahí va el código. Hoy nos enfocaremos en R, pero esta presentación se hizo en RMarkdown. Pueden correr líneas de su secuencia de códigos al usar Command+Return en Mac, o Control+Enter en la PC. También hay un botón para correr líneas o secciones enteras, así como el script entero.
1 + 2[1] 313 / 2[1] 6.52 ^ 6[1] 645 * (2 + 3)[1] 25sqrt(81)[1] 93.1 Asignando valores a objetos
Claro, R es mucho más que una calculadora básica. La computación con R incluye asignar valores a objetos, y esto se puede hacer con dos maneras básicamente equivalentes:
x = 4
x <- 4En ambos casos, x representará 4 y R guardará ese significado para las líneas subsiguientes, a menos que reasignes el valor de x.
x[1] 4x + 2[1] 6x = 8
x[1] 8Vale señalar que no se debe confundir la asignación de valor a un objeto (o variable) como una igualdad. Al pensar lo que vimos, debiéramos pensar que asignamos 4 a x o x recibe 4 o x tiene 4, mas no x es igual a 4. Aunque = es consistente con otros lenguajes de programación, muchos preferimos <- al hacer la acción tomada a cabo más evidente. Si tenían curiosidad, se prueba la igualdad con doble signo de igualdad (==), y eso produce algo distinto:
2 == 2[1] TRUE2 == 3[1] FALSEEstá bien usar nombres de variables como x para ejemplos matemáticos simples como los anteriores. Sin embargo, cuando escribas código para realizar análisis, debes tener cuidado de usar nombres descriptivos. El código donde las cosas se llaman id_sujeto, condición o edad serán más largos que x, y y z, pero tendrán mucho más sentido al volver a ellos meses después al hacer sus artículos o trabajos de investigación. Hay, eso sí, ciertas reglas para nombres: pueden ser cualquier carácter alfanumérico, pero el primer carácter ha de ser una letra. No se permiten espacios: la computadora no entiende que quieres trabajar con dos palabras como si fuera un concepto; son dos términos separados. En ese caso considera unir palabras con _ y .. R también puede trabajar con datos categóricos, no sólo numéricos:
y<-"Puerto Rico"
y[1] "Puerto Rico"Noten las comillas. ¿Qué pasa si no llevara comillas?
También podemos asignar valores lógicos como cierto, TRUE y falso, FALSE:
alive <- TRUE
asleep <- FALSEPueden comparar también valores numéricos con >, <, !=, <= y >=, que devolverán valores lógicos TRUE o FALSE.
2 < 3[1] TRUE3 <= 3[1] TRUE3 != 4 #-> noten: el símbolo complejo "!=" significa "no es igual a".[1] TRUEUna vez hayan asignado valor numérico a objetos o variables, podrán hacer cálculos:
psic <- 2
soci <- 1
econ <- 2
cipo <- 2
antr <- 1
geog <- 2
otro_cs <- 1
otras_f <- 0
taller_n = psic + soci + econ + cipo + antr + geog + otro_cs + otras_f
print(paste0("La cantidad de participantes en el taller hoy es ",taller_n, ".")) # cantidad de participantes del taller[1] "La cantidad de participantes en el taller hoy es 11."3.2 Funciones
Las funciones en R son útiles para operaciones complejas. Toman en sí insumos de argumentos o parámetros, hacen su operación y devuelven unas salidas o resultados. Ustedes llaman a la función al escribir el nombre seguido de paréntesis, con los argumentos necesarios. Por ejemplo print() y paste0 arriba. También podemos crear funciones propias:
mediana_mín_máx <- function(x){
qs <- quantile(x, c(0.5, 0, 1))
data.frame(mediana = qs[1], mín = qs[2], máx = qs[3])
}Hay varias funciones en R básico que son útiles en matemáticas:
abs(-4)[1] 4sqrt(64)[1] 8log(1.75)[1] 0.5596158Normalmente usaremos c(), la función de concatenación. Esta toma una secuencia de argumentos y la encadena en un vector. La mayoría de las funciones de estadísticas descriptivas esperan recibir al menos un vector, o algo similar a ello.
a <- c(2, 5, 7)
print(a)[1] 2 5 7#mediana_mín_máx(a)
cat(a)#, fill = T) #concatena e imprime, menos complejo que print al no dar line feeds a menos que se explicite fill=T como argumento.2 5 7sum(a) #suma[1] 14mean(a) #media[1] 4.666667sd(a) #desviación estándar[1] 2.516611var(a)[1] 6.3333334 Importando datos
Al recolectar datos en nuestros estudios, o al recibir datos pre-existentes de archivos estadísticos, es probable que nos encontremos con un archivo .csv, donde cada fila tenga la respuesta de un participante o un país, y una columna represente una variable, concepto o pregunta. Queremos importar esto a R para manipular estos datos (renombrarlos, crear nuevas variables de las existentes, fusionar conjuntos de datos que tengan puntos en común) y generar posiblemente estadísticas que resuman la información, así como analizar y visualizar las tendencias halladas.
#install.packages("tidyverse")
library(tidyverse)── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ dplyr 1.1.4 ✔ readr 2.1.5
✔ forcats 1.0.1 ✔ stringr 1.5.2
✔ ggplot2 4.0.0 ✔ tibble 3.3.0
✔ lubridate 1.9.4 ✔ tidyr 1.3.1
✔ purrr 1.1.0
── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag() masks stats::lag()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors# escribiendo un csv
#write_csv(mtcars, "mtcars.csv")
# load a csv file
d <- read_csv("mtcars.csv")Rows: 32 Columns: 11
── Column specification ────────────────────────────────────────────────────────
Delimiter: ","
dbl (11): mpg, cyl, disp, hp, drat, wt, qsec, vs, am, gear, carb
ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.Esta importación es sencilla también con Stata o SPSS usando el paquete haven:
library(haven)
haven::write_dta(mtcars, "mtcars.dta")
d2 <- haven::read_dta("mtcars.dta")
haven::write_sav(mtcars, "mtcars.sav")
d3 <- haven::read_sav("mtcars.sav")Una vez hayan importado los datos es probable que quieran echarle un vistazo a los datos, asegurarse que todo entró adecuadamente:
#mira los nombres de las variables
names(d) [1] "mpg" "cyl" "disp" "hp" "drat" "wt" "qsec" "vs" "am" "gear"
[11] "carb"#recoge los datos básicos de las variables en el conjunto (e.g. observaciones, datos ausentes, mínimo, máximo, media)
summary(d) mpg cyl disp hp
Min. :10.40 Min. :4.000 Min. : 71.1 Min. : 52.0
1st Qu.:15.43 1st Qu.:4.000 1st Qu.:120.8 1st Qu.: 96.5
Median :19.20 Median :6.000 Median :196.3 Median :123.0
Mean :20.09 Mean :6.188 Mean :230.7 Mean :146.7
3rd Qu.:22.80 3rd Qu.:8.000 3rd Qu.:326.0 3rd Qu.:180.0
Max. :33.90 Max. :8.000 Max. :472.0 Max. :335.0
drat wt qsec vs
Min. :2.760 Min. :1.513 Min. :14.50 Min. :0.0000
1st Qu.:3.080 1st Qu.:2.581 1st Qu.:16.89 1st Qu.:0.0000
Median :3.695 Median :3.325 Median :17.71 Median :0.0000
Mean :3.597 Mean :3.217 Mean :17.85 Mean :0.4375
3rd Qu.:3.920 3rd Qu.:3.610 3rd Qu.:18.90 3rd Qu.:1.0000
Max. :4.930 Max. :5.424 Max. :22.90 Max. :1.0000
am gear carb
Min. :0.0000 Min. :3.000 Min. :1.000
1st Qu.:0.0000 1st Qu.:3.000 1st Qu.:2.000
Median :0.0000 Median :4.000 Median :2.000
Mean :0.4062 Mean :3.688 Mean :2.812
3rd Qu.:1.0000 3rd Qu.:4.000 3rd Qu.:4.000
Max. :1.0000 Max. :5.000 Max. :8.000 #las primeras filas
head(d)# A tibble: 6 × 11
mpg cyl disp hp drat wt qsec vs am gear carb
<dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 21 6 160 110 3.9 2.62 16.5 0 1 4 4
2 21 6 160 110 3.9 2.88 17.0 0 1 4 4
3 22.8 4 108 93 3.85 2.32 18.6 1 1 4 1
4 21.4 6 258 110 3.08 3.22 19.4 1 0 3 1
5 18.7 8 360 175 3.15 3.44 17.0 0 0 3 2
6 18.1 6 225 105 2.76 3.46 20.2 1 0 3 1#las últimas filas
tail(d)# A tibble: 6 × 11
mpg cyl disp hp drat wt qsec vs am gear carb
<dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 26 4 120. 91 4.43 2.14 16.7 0 1 5 2
2 30.4 4 95.1 113 3.77 1.51 16.9 1 1 5 2
3 15.8 8 351 264 4.22 3.17 14.5 0 1 5 4
4 19.7 6 145 175 3.62 2.77 15.5 0 1 5 6
5 15 8 301 335 3.54 3.57 14.6 0 1 5 8
6 21.4 4 121 109 4.11 2.78 18.6 1 1 4 2#los tipos de variable
str(d)spc_tbl_ [32 × 11] (S3: spec_tbl_df/tbl_df/tbl/data.frame)
$ mpg : num [1:32] 21 21 22.8 21.4 18.7 18.1 14.3 24.4 22.8 19.2 ...
$ cyl : num [1:32] 6 6 4 6 8 6 8 4 4 6 ...
$ disp: num [1:32] 160 160 108 258 360 ...
$ hp : num [1:32] 110 110 93 110 175 105 245 62 95 123 ...
$ drat: num [1:32] 3.9 3.9 3.85 3.08 3.15 2.76 3.21 3.69 3.92 3.92 ...
$ wt : num [1:32] 2.62 2.88 2.32 3.21 3.44 ...
$ qsec: num [1:32] 16.5 17 18.6 19.4 17 ...
$ vs : num [1:32] 0 0 1 1 0 1 0 1 1 1 ...
$ am : num [1:32] 1 1 1 0 0 0 0 0 0 0 ...
$ gear: num [1:32] 4 4 4 3 3 3 3 4 4 4 ...
$ carb: num [1:32] 4 4 1 1 2 1 4 2 2 4 ...
- attr(*, "spec")=
.. cols(
.. mpg = col_double(),
.. cyl = col_double(),
.. disp = col_double(),
.. hp = col_double(),
.. drat = col_double(),
.. wt = col_double(),
.. qsec = col_double(),
.. vs = col_double(),
.. am = col_double(),
.. gear = col_double(),
.. carb = col_double()
.. )
- attr(*, "problems")=<externalptr> Recuerden que para datos pre-existentes en R o paquetes, pueden buscar más información sobre los datos al poner ?[nombredatos] en la consola. mtcars es uno de estos, así que podríamos ir a verificar en la pestaña de ayuda sobre los datos y el significado de estos.
help(mtcars)4.1 Tipos de datos
Hay cuatro tipos de datos en R. Sabiendo de ellos podemos entender lsa limitaciones de análisis de cada uno, y también podrán entender mejor algunos errores que podrían recibir. Estos son:
- Numéricos (números, enteros, dobles (aceptan los decimales))
- Caracterers (strings)
- Lógicos (C/F)
- Factores (niveles discretos; e.g., categorías)
Si ponemos la función str(d) notarán que cada variable tiene una asignatura de tipo. Pueden cambiar el tipo de datos, por ejemplo hacia categórico usando la función as.factor(). Ejemplo:
str(d)spc_tbl_ [32 × 11] (S3: spec_tbl_df/tbl_df/tbl/data.frame)
$ mpg : num [1:32] 21 21 22.8 21.4 18.7 18.1 14.3 24.4 22.8 19.2 ...
$ cyl : num [1:32] 6 6 4 6 8 6 8 4 4 6 ...
$ disp: num [1:32] 160 160 108 258 360 ...
$ hp : num [1:32] 110 110 93 110 175 105 245 62 95 123 ...
$ drat: num [1:32] 3.9 3.9 3.85 3.08 3.15 2.76 3.21 3.69 3.92 3.92 ...
$ wt : num [1:32] 2.62 2.88 2.32 3.21 3.44 ...
$ qsec: num [1:32] 16.5 17 18.6 19.4 17 ...
$ vs : num [1:32] 0 0 1 1 0 1 0 1 1 1 ...
$ am : num [1:32] 1 1 1 0 0 0 0 0 0 0 ...
$ gear: num [1:32] 4 4 4 3 3 3 3 4 4 4 ...
$ carb: num [1:32] 4 4 1 1 2 1 4 2 2 4 ...
- attr(*, "spec")=
.. cols(
.. mpg = col_double(),
.. cyl = col_double(),
.. disp = col_double(),
.. hp = col_double(),
.. drat = col_double(),
.. wt = col_double(),
.. qsec = col_double(),
.. vs = col_double(),
.. am = col_double(),
.. gear = col_double(),
.. carb = col_double()
.. )
- attr(*, "problems")=<externalptr> d$am <- as.factor(d$am)
str(d$am) Factor w/ 2 levels "0","1": 2 2 2 1 1 1 1 1 1 1 ...En estos datos am nos informa si el vehículo usa transmisión automática (0) o manual (1). Estas son categorías en realidad, representadas por una variable dummy, pues es mejor reclasificarlas de ser números continuos a categorías: al usar as.factor() le dijimos a R que am era un factor. Puedes verificar manualmente el estado o transformación de una variable al usar str(nombre_de_variable).
4.2 El operador de accesso $
A menudo queremos no entender todos los datos (pueden ser muchos) sino entender algunas variables en específico. Podemos usar un código que permite acceder a la variable dentro de un objeto: datos$nombre_var–el conjunto de datos, un signo de peso o dólar, y el nombre de la variable. Por ejemplo, d$cyl es decirle a R “dentro del conjunto de datos d, la variable cyl”. Es importante en este punto especificar el conjunto que queremos acceder, pues tenemos varias opciones.
4.2.1 Apegándonos de un conjunto de datos
Es posible que quieran trabajar con sólo uno o mayormente uno de estos conjuntos de datos. Ahí podríamos usar la opción de attach(datos) mientras operas con esos datos, y luego detach() al terminar con ellos. En este caso, le decimos a R que asuma que al llamar la variable, nos referimos a los datos que ‘pegamos’ a la memoria.
Si bien esto puede ser conveniente, podría llevar a ciertos errores por olvido o por cambios en los datos, debe ser usado sin mucha frecuencia.
Aquí ‘pegamos’ ese conjunto de datos y hacemos una operación con la variable de automático o manual am:
attach(d)The following object is masked from package:ggplot2:
mpgam <- as.integer(am)
str(am) int [1:32] 2 2 2 1 1 1 1 1 1 1 ...mean(am)[1] 1.40625sd(am)[1] 0.4989909range(am)[1] 1 2am <- as.factor(am)
str(am) Factor w/ 2 levels "1","2": 2 2 2 1 1 1 1 1 1 1 ...#sd(am) #dará error: no se puede aplicar esto a valores categóricosSi optaron por usar la función de attach, ¡asegúrense de ‘despegar’ el conjunto de datos al terminar las operaciones que iban a usar!
detach(d)4.2.2 Verificando las variables
Ahora que podemos acceder a las variables, podemos explorarlas con varias funciones existentes en R y paquetes adicionales. También podríamos crear nuestras propias funciones. Abajo algunas:
tapply(d$mpg, d$gear, mean) #Promedio de millas por galón (mpg) según la cantidad de cambios (3, 4 o 5) 3 4 5
16.10667 24.53333 21.38000 sapply(d, mean) #Promedio de cada variable del conjuntoWarning in mean.default(X[[i]], ...): argument is not numeric or logical:
returning NA mpg cyl disp hp drat wt qsec
20.090625 6.187500 230.721875 146.687500 3.596563 3.217250 17.848750
vs am gear carb
0.437500 NA 3.687500 2.812500 summary(d$mpg) #Estadísticas resumidas básicas Min. 1st Qu. Median Mean 3rd Qu. Max.
10.40 15.43 19.20 20.09 22.80 33.90 table(d$carb, d$am) #Tabla de frecuencias
0 1
1 3 4
2 6 4
3 3 0
4 7 3
6 0 1
8 0 1addmargins(table(d$carb, d$am, dnn=c('núm. de carburadores','transmisión'))) #Añade totales por fila y columna transmisión
núm. de carburadores 0 1 Sum
1 3 4 7
2 6 4 10
3 3 0 3
4 7 3 10
6 0 1 1
8 0 1 1
Sum 19 13 32# Visualizando la distribución
hist(d$wt, col = 'red3', xlab = NA, main = 'Distribución del peso')
boxplot(d$wt, ylab = 'Peso (en miles)')
plot(d$wt~d$mpg,xlab="Millas por galón",ylab='Peso (en miles)')
Nota que para estas funciones lo que se requiere en realidad es el nombre de la variable como argumento. ¿Qué creen signifiquen las opciones col =, xlab =, main =, y ylab =?
Podríamos también agrupar la información como ocurre en estos casos abajo del paquete psych.
#install.packages("psych")
library(psych)
Attaching package: 'psych'The following objects are masked from 'package:ggplot2':
%+%, alphadescribe(d$hp) # la función "describe" viene en el paquete de psych; da datos adicionales vars n mean sd median trimmed mad min max range skew kurtosis se
X1 1 32 146.69 68.56 123 141.19 77.1 52 335 283 0.73 -0.14 12.12describeBy(d$mpg, d$am) # del paquete "psych"
Descriptive statistics by group
group: 0
vars n mean sd median trimmed mad min max range skew kurtosis se
X1 1 19 17.15 3.83 17.3 17.12 3.11 10.4 24.4 14 0.01 -0.8 0.88
------------------------------------------------------------
group: 1
vars n mean sd median trimmed mad min max range skew kurtosis se
X1 1 13 24.39 6.17 22.8 24.38 6.67 15 33.9 18.9 0.05 -1.46 1.71Una forma más complicada pero útil de verificar estadísticas es la función xtabs(), es decir tabulación cruzada. El código pide una fórmula que contraste las variables. El código es un poco más complicado pero se puede interpretar como “la suma de mpg, para cada grupo de cyl y am”. Ya que cyl tiene tres grupos y am tiene dos grupos, xtabs() devuelve una tabla de seis celdas con las sumas de mpg para cada uno de estos grupos; así podríamos por ejemplo verificar como cambian las medias dadas ciertas categorías. Noten que hemos usado, sin attach() el nombre de las variables, ya que le dijimos a la función que se enfocara en las variables dentro del objeto definido en la opción data = conjunto_de_datos (aquí, data = d).
xtabs(mpg ~ cyl + am, data = d) am
cyl 0 1
4 68.7 224.6
6 76.5 61.7
8 180.6 30.8Si queríamos el promedio de millas por galón por categoría hacemos lo siguiente:
# Sumamos mpg por categorías de cyl y am
suma_mpg <- xtabs(mpg ~ cyl + am, data = d)
# Sacamos la cantidad de observaciones por categorías de cyl y am obviando el lado izquierdo de la fórmula
cantidad <- xtabs(~ cyl + am, data = d)
# Media de mpg por categorías de cyl y am
media_mpg <- suma_mpg / cantidad
media_mpg am
cyl 0 1
4 22.90000 28.07500
6 19.12500 20.56667
8 15.05000 15.40000En este caso hemos obtenido la media de millas por galón según la cantidad de cilindros del vehículo, y si la transmisión era manual o automática.
5 Entendiendo y manipulando datos: Tidyverse
Si bien podríamos usar esta información previa para entender los datos. Digamos que queremos saber qué carro es el mejor en estos datos en cuestión de millas por galón. Podríamos hacer:
library(ggplot2)
data("mpg")
?mpg
summary(mpg) manufacturer model displ year
Length:234 Length:234 Min. :1.600 Min. :1999
Class :character Class :character 1st Qu.:2.400 1st Qu.:1999
Mode :character Mode :character Median :3.300 Median :2004
Mean :3.472 Mean :2004
3rd Qu.:4.600 3rd Qu.:2008
Max. :7.000 Max. :2008
cyl trans drv cty
Min. :4.000 Length:234 Length:234 Min. : 9.00
1st Qu.:4.000 Class :character Class :character 1st Qu.:14.00
Median :6.000 Mode :character Mode :character Median :17.00
Mean :5.889 Mean :16.86
3rd Qu.:8.000 3rd Qu.:19.00
Max. :8.000 Max. :35.00
hwy fl class
Min. :12.00 Length:234 Length:234
1st Qu.:18.00 Class :character Class :character
Median :24.00 Mode :character Mode :character
Mean :23.44
3rd Qu.:27.00
Max. :44.00 max(mpg$hwy) #sin embargo esto no nos lleva a mucha más información[1] 44i_max<-which.max(mpg$hwy)
mpg$model[i_max][1] "jetta"El vehículo con mejor millaje por galón fue el Jetta. Igualmente si trabajamos datos de criminalidad en EEUU, podríamos tener
library(dslabs)
data(murders)
summary(murders) state abb region population
Length:51 Length:51 Northeast : 9 Min. : 563626
Class :character Class :character South :17 1st Qu.: 1696962
Mode :character Mode :character North Central:12 Median : 4339367
West :13 Mean : 6075769
3rd Qu.: 6636084
Max. :37253956
total
Min. : 2.0
1st Qu.: 24.5
Median : 97.0
Mean : 184.4
3rd Qu.: 268.0
Max. :1257.0 min(murders$total) #queremos igual aquí saber el estado que menos tuvo[1] 2i_min<-which.min(murders$total)
murders$state[i_min][1] "Vermont"Y tendríamos este resultado, con Vermont como el que menos asesinatos de arma de fuego registrara en 2010.
Sin embargo esto no es suficiente quizás para todo lo que queremos hacer y se ve demasiado laborioso o largo en relación a lo que obtenemos. En el caso de los datos de asesinatos por arma de fuego, hemos buscado el estado con menos asesinatos, y vemos que hay un rango enorme entre ése y el mayor, pero ¿estamos comparando chinas con chinas?
5.1 Modificando conjuntos de datos: creando una variable
Podríamos crear una tasa de asesinatos para este último conjunto de datos, con la información ya suministrada.
murders$tasa_100k<-murders$total/murders$population *10^5
head(murders) state abb region population total tasa_100k
1 Alabama AL South 4779736 135 2.824424
2 Alaska AK West 710231 19 2.675186
3 Arizona AZ West 6392017 232 3.629527
4 Arkansas AR South 2915918 93 3.189390
5 California CA West 37253956 1257 3.374138
6 Colorado CO West 5029196 65 1.292453#en tidyverse esto se puede hacer con la función mutate()
murders<-mutate(murders,tasa=total/population*10^5)
head(murders) state abb region population total tasa_100k tasa
1 Alabama AL South 4779736 135 2.824424 2.824424
2 Alaska AK West 710231 19 2.675186 2.675186
3 Arizona AZ West 6392017 232 3.629527 3.629527
4 Arkansas AR South 2915918 93 3.189390 3.189390
5 California CA West 37253956 1257 3.374138 3.374138
6 Colorado CO West 5029196 65 1.292453 1.292453Digamos que ahora queremos ver una lista más informativa e intuitiva de los datos; por ejemplo, qué otros carros figuran con buen millaje por galón:
filter(mpg, hwy>=35) # A tibble: 8 × 11
manufacturer model displ year cyl trans drv cty hwy fl class
<chr> <chr> <dbl> <int> <int> <chr> <chr> <int> <int> <chr> <chr>
1 honda civic 1.8 2008 4 auto(… f 25 36 r subc…
2 honda civic 1.8 2008 4 auto(… f 24 36 c subc…
3 toyota corolla 1.8 1999 4 manua… f 26 35 r comp…
4 toyota corolla 1.8 2008 4 manua… f 28 37 r comp…
5 toyota corolla 1.8 2008 4 auto(… f 26 35 r comp…
6 volkswagen jetta 1.9 1999 4 manua… f 33 44 d comp…
7 volkswagen new beetle 1.9 1999 4 manua… f 35 44 d subc…
8 volkswagen new beetle 1.9 1999 4 auto(… f 29 41 d subc…En este caso hemos seleccionado las filas para visualizar que tienen millaje por galón superior a 35: hemos creado un subconjunto.
5.2 El operador pipe %>% o |>
Digamos que realmente queremos reducir la cantidad de columnas en las que queremos enfocarnos, pues los datos de mpg tienen demasiadas. Podríamos crear otro subconjunto:
tabla_nueva_mpg<-select(mpg,model,year,cty,hwy)
filter(tabla_nueva_mpg, hwy>=35) # A tibble: 8 × 4
model year cty hwy
<chr> <int> <int> <int>
1 civic 2008 25 36
2 civic 2008 24 36
3 corolla 1999 26 35
4 corolla 2008 28 37
5 corolla 2008 26 35
6 jetta 1999 33 44
7 new beetle 1999 35 44
8 new beetle 1999 29 41Y nos quedamos con cuatro columnas dándonos información puntual si esto era lo que nos interesaba. Sin embargo, podríamos evitarnos la creación de objetos intermedios al usar la función que canaliza los datos en secuencia funcional así:
mpg %>% select(model,year,cty,hwy) |> filter(hwy>=35)# A tibble: 8 × 4
model year cty hwy
<chr> <int> <int> <int>
1 civic 2008 25 36
2 civic 2008 24 36
3 corolla 1999 26 35
4 corolla 2008 28 37
5 corolla 2008 26 35
6 jetta 1999 33 44
7 new beetle 1999 35 44
8 new beetle 1999 29 41Vale la pena señalar que el pipe funcionará bien con las funciones donde el primer argumento sean los datos de entrada, que canalizamos. Las funciones de tidyr y dplyr operan así y se acoplan al pipe.
5.3 Resumiendo datos explorativamente con Tidy
En esta sección cubriremos dos funciones de tidyverse, en dplyr, summarise y group_by. La primera, summarise, ofrece el cálculo de estadísticas de resumen con un código legible e intuitivo. El segundo agrupa y resume los datos por categorías inherente en las variables existentes. Por ejemplo, quizás quisiéramos calcular el promedio y desviación estándar de los vehículos en los datos de millaje por galón, pero queremos agruparlos por fabricante, o tipo de transmisión, o cualquier otro tipo de variable:
mpg %>%
group_by(manufacturer) |>
summarise(mpg_grupal=mean(hwy),
ds_grupal=sd(hwy)) #por fabricador# A tibble: 15 × 3
manufacturer mpg_grupal ds_grupal
<chr> <dbl> <dbl>
1 audi 26.4 2.18
2 chevrolet 21.9 5.11
3 dodge 17.9 3.57
4 ford 19.4 3.33
5 honda 32.6 2.55
6 hyundai 26.9 2.18
7 jeep 17.6 3.25
8 land rover 16.5 1.73
9 lincoln 17 1
10 mercury 18 1.15
11 nissan 24.6 5.09
12 pontiac 26.4 1.14
13 subaru 25.6 1.16
14 toyota 24.9 6.17
15 volkswagen 29.2 5.32mpg |>
group_by(trans) |>
summarise(mpg_grupal=mean(hwy),
ds_grupal=sd(hwy)) #por año de manufactura# A tibble: 10 × 3
trans mpg_grupal ds_grupal
<chr> <dbl> <dbl>
1 auto(av) 27.8 2.59
2 auto(l3) 27 4.24
3 auto(l4) 22.0 5.64
4 auto(l5) 20.7 6.04
5 auto(l6) 20 2.37
6 auto(s4) 25.7 1.15
7 auto(s5) 25.3 6.66
8 auto(s6) 25.2 3.99
9 manual(m5) 26.3 5.99
10 manual(m6) 24.2 5.75Podemos también aplicar una función creada previamente:
mpg |> group_by(manufacturer) |> summarise(mediana_mín_máx(hwy))# A tibble: 15 × 4
manufacturer mediana mín máx
<chr> <dbl> <dbl> <dbl>
1 audi 26 23 31
2 chevrolet 23 14 30
3 dodge 17 12 24
4 ford 18 15 26
5 honda 32 29 36
6 hyundai 26.5 24 31
7 jeep 18.5 12 22
8 land rover 16.5 15 18
9 lincoln 17 16 18
10 mercury 18 17 19
11 nissan 26 17 32
12 pontiac 26 25 28
13 subaru 26 23 27
14 toyota 26 15 37
15 volkswagen 29 23 44Digamos que queremos obtener las tasas medias de regiones específicas en los EEUU. En este caso entra más complicación pero en Tidyverse podemos usar tanto el pipe %>% como la %in% (el %in% es una función que busca parear la entrada de un valor mapeado en otro):
murders %>%
mutate(group = case_when(
abb %in% c("ME", "NH", "VT", "MA", "RI", "CT") ~ "Nueva Inglaterra",
abb %in% c("WA", "OR", "CA") ~ "Costa del Pacífico",
region == "South" ~ "el Sur",
TRUE ~ "Otras regiones")) %>%
group_by(group) %>%
summarise(tasa_100k = sum(total)/ sum(population) * 10^5)# A tibble: 4 × 2
group tasa_100k
<chr> <dbl>
1 Costa del Pacífico 2.90
2 Nueva Inglaterra 1.72
3 Otras regiones 2.71
4 el Sur 3.636 Visualización de datos: porqués
Ver números y cadenas de caracteres que forman un conjunto de datos puede ser interesante, o no, pero normalmente no tiene tanta utilidad. Por ejemplo:
library(wooldridge)
data(wage1)
head(wage1) wage educ exper tenure nonwhite female married numdep smsa northcen south
1 3.10 11 2 0 0 1 0 2 1 0 0
2 3.24 12 22 2 0 1 1 3 1 0 0
3 3.00 11 2 0 0 0 0 2 0 0 0
4 6.00 8 44 28 0 0 1 0 1 0 0
5 5.30 12 7 2 0 0 1 1 0 0 0
6 8.75 16 9 8 0 0 1 0 1 0 0
west construc ndurman trcommpu trade services profserv profocc clerocc
1 1 0 0 0 0 0 0 0 0
2 1 0 0 0 0 1 0 0 0
3 1 0 0 0 1 0 0 0 0
4 1 0 0 0 0 0 0 0 1
5 1 0 0 0 0 0 0 0 0
6 1 0 0 0 0 0 1 1 0
servocc lwage expersq tenursq
1 0 1.131402 4 0
2 1 1.175573 484 4
3 0 1.098612 4 0
4 0 1.791759 1936 784
5 0 1.667707 49 4
6 0 2.169054 81 64¿Qué aprendemos de ver estos datos así? ¿Podemos rápidamente determinar a si años de educación se traducen a mayores ingresos? ¿Podemos determinar si afecta en algo la relación marital? Para muchos humanos, es difícil extraer información con meramente mirar a números sin contexto adicional. Pero podríamos ver algo en este gráfico
Lo mismo podríamos hacer con los datos que trabajamos la semana pasada de homicidios con armas de fuego en EEUU:
library(dslabs)
library(tidyverse)
data("murders")
tail(murders) state abb region population total
46 Vermont VT Northeast 625741 2
47 Virginia VA South 8001024 250
48 Washington WA West 6724540 93
49 West Virginia WV South 1852994 27
50 Wisconsin WI North Central 5686986 97
51 Wyoming WY West 563626 5No podemos determinar con facilidad a qué estado le toca la población más grande o pequeña, y si existe alguna relación entre el tamaño de población y el total de asesinatos, o de cómo varían las tasas de asesinatos por regiones de estados en la federación estadounidense. Sin embargo, eso puede responderse sin muchas palabras con el próximo gráfico.
Una imagen vale más que mil palabras dice el dicho. Sin adentrarnos a la inferencia estadística (que cubriremos en la tercera semana de esta secuencia de talleres), hemos podido comunicar relaciones y hallazgos en datos. A veces puede ser este ejercicio uno tan convincente que no requiera análisis subsiguientes.
Vivimos en una era de creciente disponibilidad de conjuntos de datos informativos y de herramientas de software, con lo cual el uso de visualizaciones ha aumentado en diversos espacios: académicos, gubernamentales, organizaciones sociales, prensa, e industrias varias.
En R, una de las principales maneras con la que trabajaremos estos análisis visuales es de la mano de ggplot2, así como con otros paquetes que ayudan a procesar esta información. Si bien existen otros métodos de graficar en R y otros programas, la preferencia de este taller es utilizar el sistema de procesamiento que ofrece tidyverse, con ggplot2 incluido.
Esta es una secuencia de tres semanas, y en esta lección continuamos con lo aprendido la semana anterior, donde terminamos con visualizaciones sencillas y con el uso de tidyverse para manejar datos. Tenemos la meta hoy de que al culminar las segundas dos horas de esta secuencia podamos:
Entender cómo usar la gramática de gráficas
Entender cómo usar
ggplot2, continuando con lo aprendido detidyversede la semana pasada.Entender cómo utilizar en varias maneras
tidycensuspara generar mapas.
En la próxima sesión, del 10 de octubre, entraremos más en wrangling de datos, así como en la inferencia estadística, con introducciones en R sobre modelos lineales, jerárquicos y longitudinales, así como sus diagnósticos posteriores.
7 ggplot2
R ofrece varias opciones para graficar, siendo útiles las capacidades incluidas en su instalación básica. Además, existen paquetes como grid, que permite un control preciso de los elementos gráficos, y lattice, que facilita la creación de gráficos multivariados y en paneles o facetas. Sin embargo, en este taller (y en los libros de referencia usados y descritos arriba) se ha optado por usar ggplot2, ya que permite a los principiantes crear gráficos complejos y estéticos mediante una sintaxis intuitiva y fácil de recordar, dividiendo los gráficos en componentes básicos.
ggplot2 destaca por su uso de una gramática de gráficos – de donde vienen las primeras dos g – que simplifica el proceso de creación de gráficos. Al aprender unos pocos componentes esenciales de esta gramática, los usuarios pueden generar una amplia variedad de gráficos con facilidad. Además, su comportamiento por defecto está diseñado para producir resultados agradables y funcionales con código conciso y legible, lo que facilita su uso por principiantes. Como factor limitante está el que está diseñado para trabajar con tablas de formato tidy (con filas con observaciones y columnas conteniendo variables), pero un número sustancial de conjuntos de datos se trabajan en ese formato, o pueden convertirse a ese formato.
7.1 Componentes de un gráfico
Hoy construiremos varios tipos de gráficos como los que vimos arriba, así como mapas informativos. Antes que todo eso, vale señalar que los gráficos se dividen en tres componentes principales. Usaré otro gráfico que creara en 2020 durante la pandemia del Coronavirus para ejemplificar estos componentes.
- Datos:
- Estoy pasando al gráfico un conjunto de datos que corté sobre fatalidad de casos de COVID-19 hasta el 31 de julio de 2020 (los datos continúan hasta 2023).
- Este es el componente de datos del gráfico.
- Geometrías:
- El gráfico es uno de líneas, útil para varias series de tiempo. Este componente es una geometría. Otras geometrías posibles son gráficos de dispersión, histogramas, diagramas de barras, densidades suavizadas, y diagrama de cajas, entre otros.
- Mapeo estético:
- El gráfico usa señales visuales para representar en el lienzo vacío con capas distintos tipos de información provista en el conjunto de datos:
- Posiciones en el eje de x (tiempo)
- Posiciones en el eje de y (tasas de mortandad observadas)
- Color (asignado por país)
- Cada línea representa la información de un país para una serie de fechas. Estos se aclaran con una etiqueta para aclararnos esa relación de línea-color-país.
- El mapeo estético depende de la geometría utilizada.
- El gráfico usa señales visuales para representar en el lienzo vacío con capas distintos tipos de información provista en el conjunto de datos:
- Observaciones adicionales:
- Ejes x e y definidos por el rango de los datos y ambos en escalas logarítmicas.
- El gráfico incluye etiquetas, un título, etiquetas de variables, nota al calce, y el tema utilizado es uno solarizado que pareciera no muy distante al utilizado por el periódico especializado “The Financial Times”.
- Volveremos luego a estos datos para construir esta imagen si nos da tiempo en el taller, y si no, tendrán disponible el cómo hacerlo para referencia.
7.2 Lienzo vacío y capas
Manteniéndonos cerca de los datos utilizados en la semana pasada, construiremos por capas la información que va en el gráfico.
murders |> ggplot()
El primer paso de crear un gráfico de ggplot es asignar los datos a un lienzo vacío. Esto no significa poblar el lienzo con esos datos, sino pasarle al programa la información inicial, como un pintor que selecciona el tema que usará para la pintura que visualiza en su mente. Esto lo hicimos al pasar el pipe (|> o |>) los datos a ggplot, y lo que ocurrió fue que al no darle más información (capas, como la pintura en un lienzo), nos quedó un recuadro gris enmarcado por un borde blanco.
En ggplot, la información entonces se suministra al lienzo por capas, y se le pueden añadir adicionales. Esto tomará la forma de código siguiente:
DATOS |> ggplot() + CAPA 1 +
CAPA 2 + … + CAPA N
CAPA 2 + … + CAPA N
Usualmente, la primera capa que añadimos define la geometría. Si queremos hacer un diagrama de dispersión, ¿qué geometría deberíamos utilizar?
Si vemos la hoja de referencia (en la carpeta del taller 2, o accesible en esta página: https://github.com/rstudio/cheatsheets/blob/main/data-visualization.pdf), vemos que la función utilizada para crear gráficos con esta geometría puntillista es geom_point.
murders |>
ggplot() +
geom_point(aes(x = population/10^6, y = total))
En este caso ya hemos dado un paso adicional y añadido una primera capa en esta obra: le indicamos que queremos una capa que tenga una geometría de puntos, y un mapeo estético que toma las coordenadas en un plano cartesiano donde el eje de equis queda definido como population/10^6, la población de los estados o Wáshington D.C., en millones; el eje de ye queda definido entonces como el total de asesinatos con armas de fuego. En el mapeo estético entonces los puntos quedan asignados a esas coordenadas. De esta manera, las distancias entre puntos, así como otras características que queramos añadir, quedan expresadas. Esto se da a través de la función aes. Esta será de las funciones que más usen al graficar.
Noten que hasta ahora hemos trabajado este lienzo sin guardarlo como objeto. Si bien esto puede funcionar bien, es posible que queramos guardar nuestro progreso y seguir añadiendo capas adicionales. En este caso, al ejecutar el comando y guardarlo en objeto, el programa no nos dará automáticamente una actualización del gráfico; tendremos que llamar al objeto para que aparezca:
p<-murders |>
ggplot() +
geom_point(aes(x = population/10^6, y = total),size=2)
p
Tenemos entonces nuestro gráfico de dispersión inicial. Quizás no es tan informativo pero vemos una serie de puntos con el mapeo estético inicial. Noten que podríamos quitar la x= y la y= y no pasaría nada, ya que en ausencia de esta especificación, el programa entiende por defecto que lo primero que se le asigna es la información del eje horizontal, y en segundo orden el vertical:
p+
geom_text(aes(population/10^6, total, label = abb))
Aquí he añadido una geometría nueva: texto. Hay dos geometrías para esto, geom_label y geom_text, uno con el texto enmarcado en un recuadro y el segundo sin ello. Ya que cada punto (cada jurisdicción que es realmente parte de los Estados Unidos de América en este caso) tiene una etiqueta, necesitamos un mapeo estético para hacer la conexión entre los puntos y las etiquetas, así que se le asignó el mismo tal que el texto cayera exactamente en la misma coordenada que el punto. Pero esto se puede corregir:
p+
geom_text(aes(population/10^6, total, label = abb),nudge_x = 1.5) 
#nudge_y desplazaría el texto en el eje de yEn este caso hemos empujado a través del eje de equis la etiqueta de texto con un valor numérico de 1.5. Valores mayores aumentarían la distancia del texto y el punto, mientras que menores harían lo contrario.
Ahora podremos seguir con una complicación en el mapeo estético: queremos añadir una capa de color al lienzo que represente regiones de estas jurisdicciones. Esto se hace al añadir la opción de colour y dándole una variable, en este caso region.
# Crear el gráfico base con puntos coloreados por región
p <- murders |>
ggplot() +
geom_point(aes(x = population/10^6, y = total, colour = region), size = 2)
# Mostrar el gráfico
p
Ahora paso nuevamente a añadir etiquetas abreviadas, con el color por región.
# Añadir etiquetas desplazadas en el eje x (con color por región)
p<-p +
geom_text(aes(x = population/10^6, y = total, label = abb, colour = region), nudge_x = 1.5)
p
Y para darle más información al gráfico le añadimos etiquetas y títulos que sean más informativas que los nombres de variables.
#démosle más información a la gráfica de dispersión.
p+labs(
x = "Población (en millones)", # Cambia el nombre del eje x
y = "Homicidios con arma de fuego", # Cambia el nombre del eje y
colour = "Región", # Cambia el nombre de la leyenda de color
title = "Homicidios con arma de fuego en EEUU, 2010" # Título del gráfico
)
Hasta ahora hemos ido añadiendo capas pero notamos que tenemos una gran concentración de puntos en la parte inferior izquierda del lienzo: la mayoría de las jurisdicciones tienen menos de 10 millones de habitantes y menos de 400 homicidios. Esto hace leer e interpretar lo que sucede en para estos casos difícil. Podríamos entonces representar el gráfico con una transformación logarítmica al re-escalar con scale_x_continuous y scale_y_continuous:
p2<-murders |>
ggplot() +
geom_point(aes(x = population/10^6, y = total, colour = region), size = 3)
p2<- p2 + geom_text(aes(x = population/10^6, y = total, label = abb),nudge_x = 0.05)
p2<-p2+scale_x_continuous(trans = "log10")+
scale_y_continuous(transform = "log10")
p2<-p2+
labs(
x = "Población (millones, escala log.)", # Cambia el nombre del eje x
y = "Homicidios por arma de fuego (escala log.)", # Cambia el nombre del eje y
colour = "Región", # Cambia el nombre de la leyenda de color
title = "Homicidios por arma de fuego vs población, por región" # Título del gráfico
)
p2
Podríamos añadir una capa temática, usando el paquete ggthemes, que tiene un catálogo estético variado, que recomiendo verifiquen a través de https://yutannihilation.github.io/allYourFigureAreBelongToUs/ggthemes/. Esto es en adición a los que ya ofrece el paquete de ggplot2 https://ggplot2.tidyverse.org/reference/ggtheme.html.
En este caso le añadiré una visualización al estilo del semanario británico The Economist.
library(ggthemes)
p2 + theme_economist()
7.3 El gráfico que buscábamos
Normalmente queremos añadir formas o anotaciones a las figuras que no se derivan directamente del mapeo estético; algunos ejemplos incluyen etiquetas, cuadros, áreas sombreadas y líneas. Si queremos añadir una línea que represente la tasa promedio de asesinatos en todos los Estados Unidos, tendremos que determinarlo aparte con la ayuda de dplyr (parte de tidyverse), y tendremos que hacer la transformación adecuada también (logarítmica):
library(ggthemes)
library(ggrepel)
library(dslabs)
t <- murders |>
summarise(tasa = sum(total) / sum(population) * 10^6) |>
pull(tasa)
murders |>
mutate(región=case_when(
region == "Northeast" ~ "Noreste",
region == "North Central" ~ "el Midwest",
region == "West" ~ "Oeste",
region == "South" ~ "el Sur"))|>
ggplot(aes(population/10^6, total)) +
geom_abline(intercept = log10(t), lty = 2, color = "darkgrey") +
geom_point(aes(col = región), size = 3) +
geom_text_repel(aes(label = abb)) +
scale_x_log10() +
scale_y_log10() +
labs(title = "Homicidios en EEUU con arma de fuego en 2010",
x = "Población (en millones, escala log.)",
y = "Homicidios por arma de fuego (en escala log.)",
color = "Región") +
theme_linedraw()
He aquí los pasos entonces que necesitábamos para recrear la imagen arriba.’
7.4 Gráficos rápidos: qplot()
Si bien esto ha sido útil y han podido ver cómo generar gráficas complejas siguiendo la gramática de gráficos, es probable que en algún momento quieran ver rápidamente unas relaciones visuales de manera instantánea, con la intención luego de volver y darle mejor forma y seguir las complicaciones de un gráfico de ggplot2 con todas sus capas. La opción que da ggplot2 a esto es el uso de qplot().
Por ejemplo, aquí modifiqué los grupos de jurisdicciones de Estados Unidos a unos grupos específicos, y creo un gráfico para ver densidades con relleno (fill=grupo) y con líneas que también sirven para demarcar estos grupos de manera distinta (e.g. línea sólida, línea entrecortada). Noten que ha corrido todo el gráfico de la última línea.
murdersg<-murders |>
mutate(tasa= total/population *10^5,
,grupo = case_when(
abb %in% c("ME", "NH", "VT", "MA", "RI", "CT") ~ "Nueva Inglaterra",
abb %in% c("WA", "OR", "CA") ~ "Costa del Pacífico",
region == "South" ~ "el Sur",
TRUE ~ "Otras regiones"))
qplot(tasa, data = murdersg, geom= "density", fill = grupo, linetype=grupo)Warning: `qplot()` was deprecated in ggplot2 3.4.0.
En las siguientes líneas pido y ejecuto variaciones de gráficos con qplot con los datos de altura del libro de Irizarry:
data(heights)
b<-heights|>ggplot()
qplot(sex, height, data = heights, geom= "boxplot", fill = sex)
qplot(sex, height, data = heights, geom= "violin", fill = sex)
qplot(sex, height, data = heights, geom = "dotplot",
stackdir = "center", binaxis = "y", dotsize = 0.3)Bin width defaults to 1/30 of the range of the data. Pick better value with
`binwidth`.
qplot(height, data = heights, geom = "density", fill = sex)
qplot(height, data = heights, geom = "density", color = sex, linetype = sex)
Nuevamente podríamos mostrar información adicional, como medias añadidas fuera de las estéticas definidas dentro del mapeo:
mu_alt<-heights |>
group_by(sex) |>
summarise(media=mean(height))
c<-heights|>ggplot(aes(x = height))
c+geom_density(aes(fill = sex), alpha=0.4) #el alpha le da un nivel de transparencia 
#¿qué creen que pase si suben el valor a 0.9 o lo reducen a 0.1?Noten que si tenemos dos categorías, al pedir color, tiende a establecer a uno con un rojo magenta y al segundo con un azul cian (aciano o ciano), en este caso suavizado y transparentado para poder ver las distribuciones. Sin embargo, podemos editar el color manualmente, a la vez que añadimos complicaciones como una media de alturas por grupo.
c+ geom_density(aes(color = sex)) +
geom_vline(data=mu_alt, aes(xintercept=media, color=sex),
linetype="dashed") +
scale_color_manual(values=c("#999999", "#E69F00")) #noten que aquí le di los colores con código alfanumérico, un gris plateado y un dorado.
Hay varios métodos para asignar colores al nombrarlos en R. Se puede hacer por nombre (en inglés), con códigos alfanúmericos, y se puede también importar paletas de colores con paquetes específicos. Considero útil revisar las opciones en este enlace: https://r-graph-gallery.com/ggplot2-color.html
7.5 Alcanzando visualización de densidades y sus cambios a través del tiempo
library(dslabs)
data(gapminder)
gapminder |> as_tibble()# A tibble: 10,545 × 9
country year infant_mortality life_expectancy fertility population gdp
<fct> <int> <dbl> <dbl> <dbl> <dbl> <dbl>
1 Albania 1960 115. 62.9 6.19 1636054 NA
2 Algeria 1960 148. 47.5 7.65 11124892 1.38e10
3 Angola 1960 208 36.0 7.32 5270844 NA
4 Antigua… 1960 NA 63.0 4.43 54681 NA
5 Argenti… 1960 59.9 65.4 3.11 20619075 1.08e11
6 Armenia 1960 NA 66.9 4.55 1867396 NA
7 Aruba 1960 NA 65.7 4.82 54208 NA
8 Austral… 1960 20.3 70.9 3.45 10292328 9.67e10
9 Austria 1960 37.3 68.8 2.7 7065525 5.24e10
10 Azerbai… 1960 NA 61.3 5.57 3897889 NA
# ℹ 10,535 more rows
# ℹ 2 more variables: continent <fct>, region <fct>Vemos que estos datos incluyen una variedad de información a nivel de país, que incluye un número de años para ellos.
Podemos ver los grupos de países por continentes usando una matriz de gráficos, con cuadrículas paradas y acostadas:
filter(gapminder, year%in%c(1962, 2012)) |>
ggplot(aes(fertility, life_expectancy, colour = continent)) +
geom_point() +
facet_grid(continent~year)
filter(gapminder, year%in%c(1962, 2012)) |>
ggplot(aes(fertility, life_expectancy, col = continent)) +
geom_point() +
facet_grid(. ~ year)
years <- c(1962, 1980, 1990, 2000, 2012,2015)
continents <- c("Europe", "Asia")
gapminder |>
filter(year %in% years & continent %in% continents) |>
ggplot( aes(fertility, life_expectancy, colour = continent)) +
geom_point() +
facet_wrap(~year)+labs(x="Fertilidad", y="Esperanza de vida",colour="Continente")
En este caso podemos apreciar mejor cómo la dispersión y distribución de países entre los años 1960s y el presente reciente ha ido encaminado a una convergencia en estándares de vida, dejando atrás los puntos donde originaron estereotipos que aún persisten sobre cómo son las vidas y experiencias de distintas regiones del planeta (sin desestimar diferencias e inequidades presentes).
gapminder <- gapminder |> mutate(dólares_diarios = gdp/population/365)
antaño <- 1970
gapminder |>
filter(year == antaño & !is.na(gdp)) |>
mutate(region = reorder(region, dólares_diarios, FUN = median)) |>
ggplot(aes(dólares_diarios, region)) +
geom_point() +
scale_x_continuous(trans = "log2")
gapminder <- gapminder |>
mutate(grupo = case_when(
region %in% c("Western Europe", "Northern Europe","Southern Europe",
"Northern America",
"Australia and New Zealand") ~ "Occidente",
region %in% c("Eastern Asia", "South-Eastern Asia") ~ "Asia oriental",
region %in% c("Caribbean", "Central America",
"South America") ~ "Latinoamérica",
continent == "Africa" &
region != "Northern Africa" ~ "África subsahariana",
TRUE ~ "Otros"))
#le damos orden a los niveles
gapminder <- gapminder |>
mutate(grupo = factor(grupo, levels = c("Otros", "Latinoamérica",
"Asia oriental", "África subsahariana",
"Occidente")))
p <- gapminder |>
filter(year == antaño & !is.na(gdp)) |>
ggplot(aes(grupo, dólares_diarios)) +
geom_boxplot() +
scale_y_continuous(trans = "log2") +
xlab("") +
theme(axis.text.x = element_text(angle = 90, hjust = 1))
p
p + geom_point(alpha = 0.5) #añadiendo puntos para ver mejor distribuciones y cantidad de casos
library(ggridges)
p <- gapminder |>
filter(year == antaño & !is.na(dólares_diarios)) |>
ggplot(aes(dólares_diarios, grupo)) +
scale_x_continuous(trans = "log2")
p + geom_density_ridges()Picking joint bandwidth of 0.648
p + geom_density_ridges(jittered_points = TRUE,
position = position_points_jitter(height = 0),
point_shape = '|', point_size = 3,
point_alpha = 1, alpha = 0.7) Picking joint bandwidth of 0.648
antaño <- 1970
año_presente <- 2010
years <- c(antaño, año_presente)
gapminder |>
filter(year %in% years & !is.na(gdp)) |>
mutate(west = ifelse(grupo == "Occidente", "Occidente", "El Resto")) |>
ggplot(aes(dólares_diarios)) +
geom_histogram(binwidth = 1, color = "black") +
scale_x_continuous(trans = "log2") +
facet_grid(year ~ west)
Sabemos que muchos países surgieron después de 1970 (razón por la cual el histograma se nutre para el grupo El Resto), así que para comparar los países que tienen toda la información y han existido consistentemente entre esos periodos y ver si hay cambios y de dónde surgen en distribuciones podemos hacer lo siguiente:
country_list_1 <- gapminder |>
filter(year == antaño & !is.na(dólares_diarios)) |>
pull(country)
country_list_2 <- gapminder |>
filter(year == año_presente & !is.na(dólares_diarios)) |>
pull(country)
country_list <- intersect(country_list_1, country_list_2)
gapminder |>
filter(year %in% years & country %in% country_list) |>
ggplot(aes(grupo, dólares_diarios)) +
geom_boxplot() +
theme(axis.text.x = element_text(angle = 90, hjust = 1)) +
scale_y_continuous(trans = "log2") +
xlab("") +
facet_grid(. ~ year)
gapminder |>
filter(year %in% years & !is.na(dólares_diarios)) |>
ggplot(aes(dólares_diarios, grupo)) +
scale_x_continuous(trans = "log2") +
geom_density_ridges(adjust = 1.5) +
facet_grid(. ~ year)Warning in geom_density_ridges(adjust = 1.5): Ignoring unknown parameters:
`adjust`Picking joint bandwidth of 0.648Picking joint bandwidth of 0.726
gapminder |>
filter(year %in% years & country %in% country_list) |>
group_by(year) |>
mutate(weight = population/sum(population)*2) |>
ungroup() |>
ggplot(aes(dólares_diarios, fill = grupo)) +
scale_x_continuous(trans = "log2", limit = c(0.125, 300)) +
geom_density(alpha = 0.2, bw = 0.75, position = "stack") +
facet_grid(year ~ .)+labs(title = "Densidades apiladas: ingreso en dólares por día per cápita entre 1970 y 2010",
x = "Dólares por día",
y = "Densidad",
color = "Regiones")Ignoring unknown labels:
• colour : "Regiones"
8 Recapitulando
8.1 Qué aprendimos en este taller inicial
Hoy aprendimos varias cosas para empezar a usar R:
- Hablamos sobre cómo descargar e instalar R y RStudio;
- Funcionalidades posibles, como Markdown;
- Funciones y jerga específica de R;
- Cómo cargar datos y visualizarlos;
- Cómo buscar estadísticas básicas;
- Cómo llevar a cabo operaciones transformadoras a los datos;
- Empezamos a usar la sintaxis de tidyverse.
8.2 Continuamos el próximo martes
En el siguiente taller talleres que vienen continuaremos ahondando en operaciones estadísticas, visualización más avanzada, estadística inferencial y análisis de redes. Gracias por asistir hoy.